## Supplemental analysis to respond to referee points about spillover effects.
## Replicates Table 4 in Section 4 of the supplementary information.
##
## Kevin Quinn
## 10/3/2021
##

library(clubSandwich)


sink(file="SpilloverSupplement-output.txt")


mydata <- read.csv("../ScaleRaceSpring2017clean.csv")

mydata.sub <- mydata[mydata$condition %in% c("bench"),]

black.pos.1 <- grepl("mb", mydata.sub$X1a.id)
black.pos.2 <- grepl("mb", mydata.sub$X2a.id)
black.pos.3 <- grepl("mb", mydata.sub$X3a.id)
black.pos.4 <- grepl("mb", mydata.sub$X4a.id)
black.pos.5 <- grepl("mb", mydata.sub$X5a.id)
black.pos.6 <- grepl("mb", mydata.sub$X6a.id)
black.pos.7 <- grepl("mb", mydata.sub$X7a.id)
black.pos.8 <- grepl("mb", mydata.sub$X8a.id)
black.pos.9 <- grepl("mb", mydata.sub$X9a.id)
black.pos.10 <- grepl("mb", mydata.sub$X10a.id)
black.pos.11 <- grepl("mb", mydata.sub$X11a.id)
black.pos.12 <- grepl("mb", mydata.sub$X12a.id)
black.pos.13 <- grepl("mb", mydata.sub$X13a.id)
black.pos.14 <- grepl("mb", mydata.sub$X14a.id)
black.pos.15 <- grepl("mb", mydata.sub$X15a.id)
black.pos.16 <- grepl("mb", mydata.sub$X16a.id)
black.pos.17 <- grepl("mb", mydata.sub$X17a.id)
black.pos.18 <- grepl("mb", mydata.sub$X18a.id)
black.pos.19 <- grepl("mb", mydata.sub$X19a.id)
black.pos.20 <- grepl("mb", mydata.sub$X20a.id)
black.pos.21 <- grepl("mb", mydata.sub$X21a.id)
black.pos.22 <- grepl("mb", mydata.sub$X22a.id)

mydata.sub$number.black.photos.seen <- black.pos.1 + black.pos.2 +
    black.pos.3 + black.pos.4 + black.pos.5 + black.pos.6 + black.pos.7 +
    black.pos.8 + black.pos.9 + black.pos.10 + black.pos.11 + black.pos.12 +
    black.pos.13 + black.pos.14 + black.pos.15 + black.pos.16 + black.pos.17 +
    black.pos.18 + black.pos.19 + black.pos.20 + black.pos.21 + black.pos.22


target.inds <- c(23, 24)
keep.vars <- c("number.black.photos.seen",
               paste("X", target.inds, "_Q339", sep=""),
               paste("X", target.inds, "_Q329", sep=""))

mydata.sub <- mydata.sub[, keep.vars]

n.b.p.s <- NULL
resp.id <- NULL
MM <- NULL
for (i in 1:nrow(mydata.sub)){
    n.b.p.s <- c(n.b.p.s, rep(mydata.sub$number.black.photos.seen[i], 4))
    resp.id <- c(resp.id, rep(i, 4))
    MM <- c(MM, mydata.sub[i, 2], mydata.sub[i, 3], mydata.sub[i, 4],
            mydata.sub[i, 5])
}

mydata.sub.long <- data.frame(n.b.p.s=n.b.p.s, resp.id=resp.id, MM=MM)


## keep only obs with non-missing values for MM
mydata.sub.long <- na.omit(mydata.sub.long)


lm.supp.out <- lm(MM ~ n.b.p.s, data=mydata.sub.long)

tab.supp <- coef_test(lm.supp.out, vcov="CR2",
                      cluster=mydata.sub.long$resp.id,
                      test="Satterthwaite")


cat("\n")
cat("Positions 23 and 24 (just benchmark conditions)\n")
print(tab.supp)














sink()
